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Abstract 

I Backgrounds: By a new concept called "phase diagram", we compare two commonly used genotype- 
based tests for case-control genetic analysis, one is a Cochran- Armitage trend test (CAT test at x = 0.5, or 
CAT0.5) and another (called MAX2) is the maximization of two chi-square test results: one from the two-by- 
two genotype count table that combines the baseline homozygotes and hcterozygotes, and another from the 

] table that combines hcterozygotes with risk homozygotes. CAT0.5 is more suitable for multiplicative disease 
models and MAX2 is better for dominant/recessive models. Methods: Wc define the CAT0.5-MAX2 phase 

■ diagram on the disease model space such that regions where MAX2 is more powerful than CAT0.5 are separated 
from regions where the CAT0.5 is more powerful, and the task is to choose the appropriate parameterization 
to make the separation possible. Results: We find that using the difference of allele frequencies {Sp) and 
the difference of Hardy- Weinberg disequilibrium coefficients (6^) can separate the two phases well, and the 
phase boundaries are determined by the angle tan~^{Sp/Se), which is an improvement over the disease model 
selection using 5^ only. Conclusions: We argue that phase diagrams similar to the one for CAT0.5-MAX2 
have graphical appeals in understanding power performance of various tests, clarifying simulation schemes, 
summarizing case-control datasets, and guessing the possible mode of inheritance. 
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Introduction 

Comparing allele and genotype frequencies of a single marker between the patients and 
normal people (Sasieni, 1997; Lewis, 2002; Li, 2008) remains the core of a case-control genetic 
association analysis, prior to a haplotype analysis and bioinformatic analysis to determine 
the spatial extent and gene context of the signal. The success of a genetic association study 
crucially depends on study design (Amos, 2007), but the choice of test is also somewhat 
important. If we focus on genotype-based tests, i.e., each observed genotype is a sample instead 
of two allele samples (Sasieni, 1997), there are many possible types of tests to choose from. 
For example, the Pearson's goodness-of-fit test on the 2-by-3 genotype count table, Cochran- 
Armitage trend (CAT) tests with a single parameter value x which determines the relative risk 
of heterozygote with respect to the two homozygotes (Sasieni, 1997; Devlin and Roeder, 1999; 
Slager and Schaid, 2001), tests that are maximized over two or more tests (Yamada et al., 2001; 
Freidlin et al., 2002; Tokuhiro et al., 2003; Zheng et al., 2006a), and those that stepwisely use 
a few tests sequentially (Zheng et al., 2008). This quote from (Balding, 2006) might capture 
the feeling about the current state: "there is no generally accepted answer to the question of 
which single- (marker) test to use" in genetic case-control association analysis. The widespread 
use of whole-genome association in the study of complex diseases nowadays (Wellcome Trust 
Case Control Consortium, 2007) adds an urgency in having a better understanding of this 
issue. 

If one statistical test is always more powerful than another test, the first test should of course 
be used all the time. However, each test may have its own strength and weakness depending 
on the true underlying model. For example, Cochran- Armitage trend as a family of tests is to 
test the equality of genotype frequencies in the case and the control group, but as a test at a 
parameter x is to test the null hypothesis Pcasei^^) + xPcasei^CL) = Pcontroii^^) + xPcontroii^ci) 
[AA is the risk genotype, and P{AA) is its frequency; similarly for heterozygous genotype Aa). 
The CAT test at x = is most powerful for detecting recessive disease genes, whereas it is 
more powerful at x = 0.5 when Hardy- Weinberg equilibrium (HWE) holds for case and control 
groups. For a pair of tests, in the space of all possible disease models, the first test can be more 
powerful in some regions whereas the second test is more powerful in other regions. This simple 
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picture is reminiscent of the phase diagram in physics (e.g., Gibbs, 1873; Lifshitz and Landau, 
1980) where a phase can be sohd, hquid or gas depending on the temperature, pressure, or 
other relevant quantities. In our case-control genetic test example, the aim of phase diagram 
is to graphically depict regions where specific test is more powerful than others. Note that one 
should not confuse the usage of "phase" here with the meaning often used in genetics, i.e., of 
the "parental origin" of an allele in a genotype (e.g. Scheet and Stephens, 2006). 

Although the principle of our phase diagram can be established as discussed above, there 
are several practical considerations. First, a disease model is specified by many parameters, 
and one may wonder whether the phase structure can be seen in a two-dimensional space. 
Second, if the two phases are highly intermingled in one parameterization, we may ask if other 
parameterization schemes are better suited for separating the two phases. Third, besides the 
definition of phases which are based on statistical power, we can also define phases that are 
based on p-value, and the question is whether the two ways of defining phase diagrams are 
similar. This paper uses a specific pair of tests to address these questions and to show the 
utility of this approach. 

The motivation of our study is to discuss the issue of in what circumstances, one should 
use a particular test instead of others. We first study the parameterization of phase diagram, 
then construct two phase diagrams, one power-based and another p-value-based, for a pair of 
commonly used tests. The phase boundary in both versions of phase diagram will be deter- 
mined. And we examine several whole-genome association datasets using the phase diagram. 
Several other applications of the phase diagram are also discussed, including comparison of 
two random simulations of disease models, and graphical display of case-control data of many 
SNPs. 

Methods 

Case-control difference of allele frequency and Hardy- Weinberg disequilibrium 
(HWD) coefficient given a disease model 

A single-locus biallelic disease model can be specified by four parameters: three penetrances 
fi =Prob (disease I Gj) {i = 0,1,2) where {Gj} are the three genotypes aa,aA, AA, and the 
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population A-allele frequency p. Alternatively, one can use the allele frequency p, two relative 
genotype risks Aj = /j//o {i = 1, 2), and disease prevalence K = /o[(l— p)^ + 2p(l— j>)Ai+p^A2] 
to specify the disease model. Denote A-allele frequency in case and control group as pi and pq, 
and HWD coefficient (Weir, 1996) in case and control group as ei and eo; all four {pi,Po, ei, eo) 
can be derived from a given disease model, under the assumption that HWE holds true for 
the general population that contains both cases and controls (Wittke-Thompson et al., 2005). 
We use the difference of A-allele frequency {6p) and the difference of HWD coefficients (5^) in 
case and control group as the parameters for our phase diagram: 

/o(p2A2 -p)Ai) p\l - /0A2) -p)(l - /oAi) 



Sp = Pi- Po 
6e = ei- eo 



K 1-K 

/oV(l -P)'(A2 - A?) /op2(l -p)2(2Ai - 1 - A2 - /oA? + /0A2 



(1) 



When Al > 1 and A2 > 1, A-allele is enriched in the case group, and 6p is always positive. On 
the other hand, if Ai < 1 and A2 < 1, A-allele is depleted in the case group and Sp is negative. 
To be consistent, we call A the risk allele when 6p > 0. Whether A-allele is a minor allele 
(p < 0.5) or a major allele {p > 0.5) does not by itself affect the sign of 6p and 6^- 

Case-control difference of allele frequency and Hardy- Weinberg disequilibrium co- 
efficient given a genotype count table 

Table 1(A) shows a case-control genotype count table {A^y} {i =1,0 for case, control, and 
j =0,1,2 for genotypes with 0, 1, 2 copies of allele A) with total A^ = A^i + A^o samples. 
The same table is parameterized in Table 1(B) using the estimated {pi} and {ej} {i =1,0 for 
case, control). The estimated A-allele frequency {6p) and the estimated difference of HWD 
coefficients (S^) in case and control group can be obtained from Table 1(A) as: 

A^n/2 + A^i2 A^oi/2 + A^o2 



= Pi- Po 



Ni No 



, . . A^i2 /AV2 + A^V A^02 ^ /A01/2 + A/02 V 

For notation simplicity, the hat () will be removed later on, and whether the model-based or 
data-based usage is applied should be clear from the context. 
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Note that switching the A and a allele, or equivalently, switching the first and the third 
column in Table 1, changes the sign of 6p, whereas the sign of 6^ is unaffected. Eq.(l2]) also 
shows another advantage of using the difference of two HWD coefficients: if the disequilibrium 
is sensitive to typing errors, its effect is minimized when the difference 6^ is used. 

Cochran-Armitage trend test and MAX2 test 

Cochran-Armitage trend (CAT) test at parameter x assigns a score of 0, x, and 1 for 
genotypes aa, aA, and AA, for log-risk relative to the baseline aa genotype (Sasieni, 1997; 
Devlin and Roeder, 1999; Slager and Schaid, 2001). Sometimes, the name of Cochran-Armitage 
trend test is used to only refer to that at parameter x = 0.5 (Balding, 2006), which we will call 
as CAT0.5. On the other hand, CAT at parameter settings of x = and x = 1 is equivalent 
to assuming recessive and dominant models, which we will refer to as CATO and CATl. 

With the genotype count in Table 1, the test statistic of CATO. 5, CATO, and CATl can be 
derived (see, e.g., Sasieni, 1997). We re-parameterize these test statistics using a new set of 
parameters including Sp and Se'- 



X\CATO) 



AiAo {2p6p + 



N (p2 + e){l-p^ -e) 



X^CATl) = Ml W-P)K-^.?_ ,3) 

^ N (2p-p2_e)(i-2p + p2 + e) 

where 6p = p, - po, p = {Ni/N)pi + (Ao/A)po, ^ = {NjN)pl + (Ao/A)p2, e = {N,/N)ei + 
(Ao/A)eo, P= (pi +po)/2. CATO. 5 usually leads to very similar result to the allele-based test. 
But since CATO. 5 is a genotype-based test, it does not have the problem of allele-based test 
for artificially doubling the sample size (Sasieni, 1997). 

The MAX2 test is defined by the maximization of the CATO and CATl test statistics: 
X^{MAX2) = max{CATO,CATl). Ahhough MAX2 has been used in a few analyses (e.g., 
Yamada et al., 2001; Tokuhiro et al., 2003), it did not have a formal name. The name of MAX2 
used here is to distinguish it from the MAX3 test {ma.x{CAT0,CATl,CAT0.5)) proposed in 
(Freidlin et al., 2002). 
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Since MAX2 involves a multiple testing whereas CAT0.5 does not, its values are calculated 
differently. For CAT0.5, the p-value is simply derived by the distribution with 1 degree of 
freedom. For MAX2, Dunn-Sidak multiple testing correction (Ury, 1976) is exact if CATO and 
CATl test statistics were independent: 

PMAX2 = Po{X\MAX2) > T) = 1 - Po{X\MAX2) < T) = 1 - P^{X^ < Tf 

= l-{l-px2f = 2px2-p\2 (4) 

where Pq is the null distribution, px2 is the value for a x^-distributed test statistic. In R code 
( http:/ / www.r-project.org/^ , the command for calculating p- value for MAX2 is l-pchisq(MAX2, 
df=l)~2. 

For dataset generated by a known disease model, CAT test statistics follow chi-square 
distributions with non-central parameters. The non-central parameters for CATO. 5, CATO 
and CATl are given in Eq.Q, only that 6p, 5^, and other parameters are determined by the 
disease model, not from the data. Alternatively, power can be determined empirically by 
simulation. 

Simulation 

For empirical power calculation, we sampled A^^ (=5000) replicates of genotypes for Ncase 
(=500) cases and Ncontroi (=500) controls, given a disease model. In FiglJl Nm (=10000) 
disease models were generated randomly. The relative genotype risks Ai and A2 are randomly 
selected from a range (e.g., (0.5-2)), and the population A-allele frequency is randomly selected 
(e.g., from (0.1-0.9)). The type I error is set at 0.05 and we have determined the test statistic 
threshold for MAX2 either by permutation or by Dunn-Sidak formula. Due to the consistency 
between the two approaches, the type I error is controlled mostly by using the Dunn-Sidak 
formula in Eq.Q. The empirical power of CATO. 5 or MAX2 at the given type I error is 
determined by the proportion of replicates that exceeds the threshold. 

For the phase diagram in FiglHunder the null model (i.e., same allele and genotype frequency 
according to the HWE), Nr (=2000) replicates of genotype were generated for Ncase (=500) 
cases and Ncontroi (=500) controls. 
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Results 

Phase diagrams based on power given the disease model 

In the S^-6p space, known types of disease models such as dominant, recessive, additive, 
multiplicative, and over-dominant models fall in different regions of the plane, after requiring 
the risk allele A to have higher frequency in cases than in controls, as can be seen from FiglH 
For example, recessive models reside in the first quadrant, additive, dominant, over-dominant 
models are in the second quadrant, and the multiplicative models sit along the ?/-axis. 

We pick CAT0.5 and MAX2 as the two tests to compare for the following reasons. First, 
allele-based test is still one of the most commonly used tests in case-control genetic analysis, 
and we would like to choose a similar genotype-based test. Second, we want to choose a test 
that is robust against disease model mis-specification. These two considerations lead to CAT0.5 
and MAX2. Figl2] shows the power-based CAT0.5-MAX2 phase diagram, where empirically 
obtained statistical power of CAT0.5 and MAX2 are compared while controlling the type I 
error, in the space parameterized by 6^ and 6p. 

The region in Figl2]where power(MAX2) > power(CAT0.5), or phase 1, covers most of the 
model space. On the other hand, region where power(CAT0.5) > power(MAX2), or phase 2, 
is limited to a narrow angle around the ?/-axis. For regions far away from the x-axis, both 
tests lead to close to 100% power (the symbol "1" is used to mark the points when both 
power(CAT0.5) and power (MAX2) are larger than 0.99). The phase boundary in the first 
(and the third) quadrant can be roughly approximated by the line 6 = tan~^{6p/6e) = 73.125° 
(137r/32). The phase boundary in the second (and the fourth) quadrant is not sharp with 
some degree of overlap between the two phases. However, the line 6 = tan~^{5p/ 5^) = 106.875° 
(197r/32) seems to provide a reasonable boundary to phase 2 points. 

The phase structure presented in Figj2] is consistent with our current knowledge that allele- 
based test and CAT0.5 are most powerful for multiplicative models. In fact, it can be shown 
that non-central parameter in the distribution for the allele-based test is strictly larger than 
that of either CATO or CATl (i.e., either one of two ways to combine heterozygote counts with 
the homozygote counts), if (1) ei = eo = and (2) Ni = Nq (Suh and Li, 2007). However, it 
was somewhat surprising that the two phases in Figl2] are not separated by vertical lines. The 
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two-phase genetic model selection method proposed in (Zheng and Ng, 2008) attempts to infer 
the underlying disease model by 6^ value alone. Here we show that the model selection could 
be more accurate if both and Sp parameters are considered. 

To have a better understanding of the phase structure, in FigjSl we plot the power of CAT0.5 
as a function of \6p\ as well as power of MAX2 as a function of radius ^y6^'+~6^. As expected, 
the power of CAT0.5 increases as the allele frequency difference is larger, because CAT0.5 is 
very close to the allele-based test which is designed to detect allele frequency differences. On 
the other hand, it was unexpected that the power of MAX2 increases with the radius, and the 
increase follows a similar way as power of CAT0.5 increases with \6p\. If we draw the two plots 
in one, the two more or less overlap. In a crude approximation, suppose power(CAT0.5) vs. 
\6p\ and power (MAX2) vs. ^5^ + 5^ follow the same functional form, then the power of MAX2 
should be larger than that of CAT0.5 over the whole plane except the y-axis, because radius is 
always larger than the y value except for the y-axis. Although this is only an approximation, 
Figj3] helps us to understand why the robust MAX2 test is more likely to be more powerful 
than the CAT0.5 test in the space of all possible models. 

Phase diagrams based on p-values given the genotype data 

Power analysis is always discussed in a theoretical context because it requires our knowing 
of the true disease model. In reality, the disease model is not supposed to be known. Towards a 
data-driven concept of phase diagram, we define the following phase diagram based on p- values. 

As can be seen from Eq.Q, the relative magnitude of CAT0.5 and MAX2 test statistics 
depends on the estimated value of four parameters pi,Po, ei, ^o, as the common factor NiNq/N 
is canceled out. If we assume that averaged quantities such as p and e do not vary dramati- 
cally with different pi,po, ei, eo values, then the relative magnitude of CAT0.5 and MAX2 test 
statistics are mainly determined by the parameters 6p and 6^. 

Because of the multiple testing in MAX2, even if a MAX2 test statistic is larger than that 
of CAT0.5, the MAX2 p- value may not be smaller than the p- value by CAT0.5. It should 
be noted that the assumption to apply Dunn-Sidak formula in Eq.(j4]), i.e., the independence 
between the two terms to be maximized, does not hold exactly. However, the correlation 
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between CATO and CATl is actually small. It was shown in (Freidlin et al., 2002; Zheng et al., 
2006a) that the correlation between CATO and CATl test statistics under the null hypothesis 
(Pi = Po)) if HWE is true, is \/pi{l — pi)/[(l + pi){2 — Pi)]- This correlation is at most 1/3 
which is reached when Pi = Po = 0.5. Simulation shows that Dunn-Sidak formula leads to an 
accurate estimation of Pmax2, despite a weak violation of the assumption. This is in a contrast 
with MAX3, where the correlation between CATO. 5 and either CATl or CATO is too strong 
to apply the Dunn-Sidak formula (Gonzalez et al., 2008; Li et al., 2008). 

FigjH shows the phase diagram based on comparing p-values of CATO. 5 and MAX2 when 
genotype count tables are randomly sampled from the null model. The MAX2 (CATO. 5) test 
is preferred over CATO. 5 (MAX2) in phase 1 (2) because it has a smaller p-value (black (grey) 
dots). When we draw the same straight lines as Fig l2] with 73.125° and 106.875°, a similar 
observation can be made that the phase separation is much better in the first (and the third) 
quadrant, and a low level of overlap occurs in the second (in the fourth) quadrant. We conclude 
that the CAT0.5-MAX2 phase diagram based on p- values is very similar to the CAT0.5-MAX2 
phase diagram based on power. Because the ranges of x and y in FigJH are smaller than those 
in FiglJl another understanding of p-value-based phase diagram is that it is the "extension" 
of power-based phase diagram towards the origin. 

Illustration of the CAT0.5-MAX2 phase digram by results of several genome-wide 
association studies 

Recently, the largest scale whole-genome association study was carried out by Wellcome 
Trust Case Control Consortium on several common complex diseases (Wellcome Trust Case 
Control Consortium, 2007). We take the genotype count for SNPs that showed the strongest, 
and for many of them, validated association signal using extra samples, to illustrate the use of 
phase diagram. The association signal in (Wellcome Trust Case Control Consortium, 2007) s 
mainly ranked by the CATO. 5 test. 

Table 2(A) lists the raw genotype count data of 11 SNPs that are associated with one of 
these diseases: ankylosing spondylitis (Wellcome Trust Case Control Consortium & Australo- 
Anglo- American Spondylitis Consortium, 2007), type 1 diabetes (Todd et al., 2007), Crohn's 
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disease (Parkes et al., 2007), and type 2 diabetes (Zeggini et al., 2007). Table 2(B) shows 
the estimated case or control minor allele frequency {pi,po) and their difference (Sp), case or 
control HWD coefficients (ei,eo) and their difference (5^), the angle 6 in the 6p — 6^ phase 
diagram with respect to either x ot y axis, and CAT test statistics at x = 0.5, 0, f . 

There are several observations made from Table 2(B). As expected, the estimated HWD co- 
efficient in case group ei is usually larger in magnitude than that in control group eo. However, 
the largest observed ei is only around 0.01. On the other hand, allele frequency difference is 
large (0.03-0.07) due to the fact that these SNPs are selected by significant CAT0.5 test. A 
consequence of the two facts is that the angle with respect to the y-axis in Table 2(B) tends to 
be small, with the exception of SNP rs27044 which is associated with ankylosing spondylitis 
(^-90° = 16.4°). 

The closeness to the y axis of these SNPs on the phase diagram should indicate, on average, 
CAT0.5 test to be better than MAX2 test. Indeed, the ^^(CATO.S) test statistics are larger 
than X^(MAX2) except for two SNPs: rs27044 and rs2542151. It is not surprising to see 
rs27044 in the exception list as it has the largest angle with respect to the |/-axis, and the 
negative sign of 6^ indicates that the disease model for rs27044 is more likely dominant. For 
rs27044, MAX2 test leads to a more significant resuh (p-value = 1.8 xlO"'^) than the CAT0.5 
test (p-value= 10~^). For rs2542151, its position in the phase diagram forms a smaller angle 
with the y-axis than the rs27044 SNP (8.5° vs. 16.4°), but the angle is still large enough such 
that MAX2 leads to a more significant test (p-value =3.7 x 10"^^) than CAT0.5 test (p-value 
= 1.9 X 10-^^). 

Discussion and conclusions 

Our choice of 6p is to capture the linear or first-order term from the disease model or data, 
and (5e to capture the nonlinear or second-order term. In quantitative genetics, there is a 
similar approach in using the additive variance and dominance variance (Fisher, 1918; 
Falconer and Mackay, 1996). When these concepts from quantitative genetics are translated to 
dichotomous traits, al = 2p(l-p)(p(/2-/l) + (l-p)(/l-/o))^ and aj = p'{l-pYU2+fo-Vi? 
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(see, e.g., Blackwelder and Elston, 1985). Using Eq.([T]), we have 

_ Pf2 + {1-P)h-K Pih - h) + q{h - fo) 

~ Kil-K) ^" K{1 - K) 

In other words, the square of our first-order parameter 5"^ is proportional to the additive 

variance a^. Similarly, using Eq.([T]) for 5^, we have 

~ i^2(l-i^)2 1/2/0 - Jl)+ _ 1/2 + /O - ^/ij- Ibj 

When the disease prevalence is low, 5^ is roughly proportional to /2/0 — /i , as versus the 
/2 + /o — 2/1 expression in aj. Actually, the HWD coefficient in the control group, which 
is small, is proportional to the difference between /2/0 — fl and /2 + /o ~ 2/i, and two are 
approximately equal when Ai, A2 are small (Zheng et al., 2006b). Crudely, the square of 5^ can 
be said to be proportional to the dominance variance aj. 

The idea to use the test that is most powerful to the underlying model sounds straightfor- 
ward, but in reality the true disease model is unknown. There have been attempts to infer 
the disease model by the HWD information. For example, (Wittke-Thompson et al., 2005) 
distinguishes HWD from different disease models and proposed its use for data fitting (note 
that the additive model defined in (Wittke-Thompson et al., 2005), A2 = 2Ai, is different from 
that defined here, /2 — /i = /i — /o oi' -^2 = 2Ai — 1). In (Zheng and Ng, 2008), the signs of ei 
and eo are used for genetic model selection: (+,— ) for recessive models, (— ,+) for dominant 
models, and (— ,— ) for multiplicative and additive models. Since the amount of HWD in con- 
trol group is usually much smaller than that in case group, the sign of 5^ = ei — eo ~ niay 
serve the purpose in selecting recessive and dominant models. All these previous works use 
HWD alone in genetic model selection, without considering a joint effect of HWD and allele 
frequency difference. 

The result discussed in this paper shows that a joint consideration of 5^ and 5p could be 
more effective, than a consideration of 5^ only, in selecting disease model. The following 
simple procedure might be reasonable: first draw a line from origin to the data-determined 
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(x, y) = [S^, Sp) position, then calculate the angle formed by this line and the y-axis. If the angle 
is smaller than 37r/32 (or 16.875°), the underlying model is more likely to be multiplicative 
and CAT0.5 is the preferred test to use. On the other hand, if the angle with respect to 
the y-axis is larger, the underlying model is more likely to be recessive (if it's in the first 
quadrant) or dominant (second quadrant), and MAX2 is preferred over CAT0.5. A caution on 
this procedure is that, unless the sample size is very large and unless the true model is away 
from the phase transition line, the disease model can still be incorrectly inferred. 

Comparing the power of two tests is always tricky because the answer depends on what 
is known about the true model. It is very much in the spirit of Bayesian statistics that the 
posterior probability distribution depends on the prior. Any claim on discovering a more 
powerful test may contain a caveat on how the alternative model is sampled. Phase diagram 
discussed here provides a framework to visually inspect distributions on the 6p-6e, for the 
simulated models, 

FigjS] shows the distribution on the phase diagram of two different ways of simulating 
genotype count tables. The first is by randomly sampling two allele frequencies Pcase and 
Pcontroi that are close to each other (difference is less than 0.1), then the genotype frequencies 
in the model are determined by the HWE. One such model is used to generate one replicate 
dataset and S^, Sp values are determined from the replicate. The second way to generate a 
random model is to randomly sample two sets of genotype frequencies, for case and control 
group, that are close to each other (difference is less than 0.1), and that model is used to 
simulate one replicate dataset. 

When we compare the empirical power of MAX2 and MAX3, MAX3 is more powerful than 
MAX2 in the first simulation, whereas MAX2 is slightly more powerful than MAX3 in the 
second simulation. The phase diagram in FigJS] easily solves the puzzle: from FigJSl it can 
be seen that simulated datasets by the first method centered around the |/-axis, as genotype 
frequencies in the model are obtained by HWE, whereas those by the second method are more 
widespread in the x-axis. If we require a graphic showing of the phase diagram like Figj5] in 
any power comparison between two tests, there would be less misunderstanding of seemingly 
conflicting empirical results. 

The phase-diagram can also be used to summarize a case-control dataset with many SNPs. 
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Fig. in] shows 2147 SNPs on chromosome 18q21 for 460 cases and 459 controls used in Genetic 
Analysis Workshop 15 (Amos et al., 2007; Wilcox et al., 2007), with 6p and 6^ calculated by 
Eq.([2]). These SNPs are filtered from the original list of 2719 SNPs by requiring no more than 
10 missing typings and minor allele frequency to be larger than 0.05. For some SNPs, alleles 
are switched to make sure 6p is positive. One striking visual impression of Figl6] is that the 
SNPs with the largest 6p values are not located on y-axis, but in the first quadrant, indicating 
that recessive model better describes the effect of these SNPs than the multiplicative model. 

In Figjni the SNP rs3745064 at position 53.305Mb exhibits the largest Sp value (=0.076) 
(Kuo et al., 2007; Tapper et al., 2007). The distribution of case samples among the three 
genotypes is 37, 180, 243, and that of control samples is 50, 223, 186, leading to pi = 0.724, 
Po = 0.648, 6i = 0.0042, 6o = —0.015. For this SNP, we expect the largest test statistic value 
to be CAT(x = 0) because it is located in the first quadrant on the phase diagram. Indeed, 
CAT(a; = 0) = 13.97, CAT(x = 0.5)=12.47, CAT(x = 1)=2.18, and MAX2 test statistic is more 
powerful than CAT{x = 0.5). 

Using the phase diagram to examine the most recent whole-genome association data shows 
that SNPs with the strongest association signal tend to be multiplicative, and CAT0.5 test is 
more powerful than MAX2 test in this situation. However, we should not exclude the possibility 
that it is due to a selection bias, as the top ranking genes were chosen by CAT0.5 test result. 
Also, if the most significant SNPs exhibit larger allele frequency differences, whereas their 6^ 
value is limited, their positions in the phase diagram is expected to be closer to the y-axis. 

Whether the result in Table 2(B) is due to selection bias or not, our approach could be 
useful in analyzing whole-genome association data, as illustrated by FiglHl Applying MAX2 
may change the rank order of some SNPs that are near the top, and consequently change 
the pool of SNPs to be studied further. If MAX2 test does improve the p-value over CAT0.5 
for a SNP, one can check whether the HWD mainly occurs in the case instead of the control 
group. The location of a SNP on the phase diagram provides a quick filtering of SNPs where 
the inclusion of non-multiplicative models may improve the association signal, such as the 
example of rs3745064 on Fig|6l Due to a high cost of typing extra samples in validating the 
associated SNPs for the second round, it is important to carry out careful analyses including 
the consideration of alternative disease models. 



Li, Suh, Yang 



14 



In conclusion, using the phase diagram is to partition the violation of null hypothesis of 
equal genotype frequency in case and control groups into two components, one 6p for allele 
frequency difference and another 5^ for difference in HWD coefficients. The relative magnitude 
of 6p and 6^ determines which test, e.g., between CAT0.5 and MAX2, is more powerful. The 
phase diagram highlights the point that no uniformly powerful test exists, and a test is only 
more powerful regionally in the model space. We believe the use of phase diagram can aid the 
design of test when some knowledge of the mode of inheritance is available, as well as inferring 
the underlying mode of inheritance from the data. 
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(A) genotype count table 

sample size aa aA AA 

case(l) Ni Nw Nn Nu 

control(O) Nq TVqo A'^oi A^02 



(B) same genotype count table parameterized by pi,po,ei, eo 
aa aA AA 

case(l) 7Vi[(l-pi)2 + ei] 7Vi[2(l - - 2ei] Ni[pl + ci] 

control(O) iVo[(l - po)^ + eo] ^o[2(l ~ po)po - 2eo] Na[pl + eo] 



Table 1: Genotype count tabic (A) and its parameterization using pi,po, ei, eo (B). 
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(A) genotype count table 



disease gene SNP chromosonie case(aa/aA/AA) control (aa/aA/AA) 



AS" 


ARTSl 




rs27044 


5 


793/553/119 


395/432/94 




AS 


LNPEP 




rs2303138 


5 


738/176/8 


1269/193/4 




TID^ 


C12orf30 


rsl7696736 


12 


1984/3115/1168 


1545/2891/1373 




TID 


ERBB3 




rs2292239 


12 


2592/2805/801 


1956/2816/946 




TID 


KIAA0350 


rsl2708716 


16 


2652/2857/779 


2834/2429/569 




TID 


PTPN2 




rs2542151 


18 


4219/1635/182 


3628/1889/220 




TID 


PTPN22 


rs2476601 


1 


4674/998/55 


3754/1580/178 




CD'^ 


IRGM 




rsl3361189 


5 


2061/705/59 


8907/1476/54 




CD 


N0D2 




rsl7221417 


16 


1505/1175/256 


747/754/245 




CD 


IL23R 




rsll805303 


1 


1385/1236/313 


655/815/276 




T2D'' 


FTO 




rs8050136 


16 


1063/1407/464 


550/987/378 










(B) position on the phase dia^ 


;ram 








disease/SNP 


Pi 


Po 


s. 


ei 


eo 




e 


61-90° 


CAT(0.5/0/l) 


AS/rs27044 


.34 


.27 


.067 


-.011 


.0083 


-.020 


106.4 


16.4 


23.9/3.0/28.6 


AS/rs2303138 


.069 


.10 


-.036 


-.0020 


-.0022 


2 X 10-4 


-89.7 


0.3 


19.4/4.0/17.9 


TlD/rsl7696736 


.49 


.43 


.050 


9 X 10^-* 


-.0028 


.0037 


85.8 


-4.2 


61.5/45.3/37.3 


TlD/rs2292239 


.41 


.36 


.056 


-.004 


.0028 


-.0069 


97.0 


7.0 


79.3/31.2/73.0 


TlD/rsl2708716 


.31 


.35 


-.045 


.004 


6 X 10-4 


.0034 


-85.7 


4.3 


55.4/21.2/50.3 


TlD/rs2542151 


.20 


.17 


.037 


-.0029 


.0027 


-.0056 


98.5 


8.5 


54.7/5.6/58.7 


TlD/rs2476601 


.18 . 


097 


.079 


.0015 


2 X 10"'' 


.0012 


89.1 


-0.9 


292.6/71.2/273.2 


CD/rsl3361189 


.076 


.15 


-.070 ^ 


-6 X 10-'' 


-3 X 10--^ 


-3 X 10-4 


-90.2 


-0.2 


261.5/65.0/238.4 


CD/rsl7221417 


.36 


.29 


.069 


.013 


.0047 


.0088 


82.8 


-7.2 


46.5/32.3/31.5 


CD/rsll805303 


.39 


.32 


.074 


.0048 


.0060 


-.0012 


90.9 


0.9 


51.7/26.3/41.8 


T2D/rs8050136 


.46 


.40 


.057 


-.0097 


-2 X 10-"^ 


-.0095 


99.5 


9.5 


31.4/12.4/29.4 



Table 2: (A) Genotype counts of SNPs that show strong (and validated) association signal with several 
common diseases, a. ankylosing spondylitis (Wellcome Trust Case Control Consortium & Australo-Anglo- 
American Spondylitis Consortium, 2007). b. type 1 diabetes (Todd ct al.. 2007). c. Crohn's disease (Parkes et 
al., 2007). d. type 2 diabetes (Zeggini et al., 2007). Other information includes gene and SNP name, as well as 
the chromosome number on which the SNP is located. (B) For the SNP listed in (A), minor allele frequency 
in case (pi) and control (po) group, their difference (5p = pi —po), HWD coefficients in case (ei) and control 
(eo) group, and their difference (S^ = ei — eo)j angle 6 ~ tan^^{Sp/S^) in the phase diagram, angle with respect 
to the y-axis {0 - 90° if 6* > 0, or 6* + 90° if 6* < 0), and test statistics for CAT0.5, CATO, CATl. 
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Figure 1 : Locations in the phase diagram parameterized by x = (5^ and y = Sp for some commonly encountered 
disease models, recessive (R/r), dominant (D/d), additive (A/a), multiplicative (M/m) and over-dominant 
(O/o). The upper case is used when the power of MAX2 is larger than the power of CATO.S, and lower case for 
the opposite. The model parameters are sampled from these ranges, population A-allclc frequency is randomly 
sampled from (0-0.5), phenocopy rate /o from (10^^-10^'^), with the exception of over-dominant models, A2 
from (1.001-2) and Ai is derived from A2 (Ai = A2, 1, (A2 + l)/2, ^/X2 for for dominant, recessive, additive 
models, and multiplicative models). For over-dominant models, Ai is randomly chosen from (1.001-2), and A2 
from (1.001, Ai). 
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Figure 2: The CAT0.5-MAX2 phase digram, with 5^ = ei — eo (difference of two HWD coefficients in case 
and control group) as the z-axis, and i5p = pi — pa (difference of allele frequencies in case and control group) 
as the y-axis, determined by comparing the empirical power of CAT0.5 and MAX2 tests. More than 10000 
disease models are randomly sampled, each is used to randomly generate 5000 replicates of 500 cases and 500 
controls. The empirical power of CAT0.5 and MAX2 based on these 5000 replicates are compared when the 
type I error is controlled at 0.05. A black circle is drawn when the empirical power of MAX2 is larger than 
that of CAT0.5, and a gray circle when CAT0.5 is more powerful than MAX2. If the empirical powers of both 
tests are larger than 0.99, a symbol "1" is marked. The two dashed lines form angles of 73.125° and 106.875° 
(137r/32 and 197r/32) with the a--axis, or, ±16.875° (37r/32) with the y-axis. 
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Figure 3: (A) The statistical power of MAX2 as a function of the radius in the 6t — 6p plane, ^S'^ + 5^ . (B) 
The statistical power of CAT0.5 as a function of vertical distance to the x-axis, \dp\. 
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Figure 4: The CAT0.5-MAX2 phase diagram based on p-value of two tests. 2000 genotype count tables are 
simulated using the null model with 500 cases and 500 controls, and A-allele frequency being samples from the 
range of (10^^-0.5). If the p-value by MAX2 is smaller (larger) than that by CAT0.5, a black (grey) circle is 
drawn. The two dashed lines are the same as those in FiglJl with angles of 73.125° and 106.875° (137r/32 and 
197r/32) with the x-axis. 
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Figure 5: Distribution of datascts in two random sampling of disease models: (1) pA,case is randomly selected 
from (0.1, 0.9); and pA.controi from (pA.case-0.1, pA,case+0.1). Genotype frequencies in cases and in controls 
follow the Hardy- Weinberg equilibrium within the respective group. One replicate of genotype counts is 
simulated which is used to calculated the location in phase diagram by Eq.(l2]). (2) Genotypes aa, Aa, AA in 
case group are randomly selected by the interval probabilities: from to Li, from Li to L2, and from L2 to 1, 
where Li and L2 are sampled from (0.1,0.45) and (0.55, 0.9). The genotype frequencies for control group are 
determined by the interval probabilities [0, L'l, L'2, 1] also, with L'l and L2 randomly sampled from (Li-0.1, 
Li+0.1), (L2-O.I, ^2+0.1). 
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Figure 6: Distribution of 2147 SNPs on chromosome 18q of a case-control dataset (Amos et al., 2007) on the 
6e-Sp phase diagram. The SNP with the largest allele frequency difference, rs3745064, is marked with a circle. 
Note that rs3745064 is not on the y-axis, but inside the first quadrant, indicating that the recessive model 
describes its effect better than the multiplicative model. 



